M O N T E Carlo Simulation of the Total Radial Distribution Function for Interlayer Water in Sodium a N D Potassium Montmorillonites

-Monte Carlo simulations based on tested water-water, cation-water, and water-clay potential functions were applied to calculate radial distribution functions for O-O, O-H and H-H spatial correlations in the interlayer region of the two-layer hydrates of Naand K-montmorillonite. The simulated radial distribution functions then were used to compute the total radial distribution function for interlayer water, a physical quantity that can be determined experimentally by H/D isotopic-difference neutron diffraction. The simulated total radial distribution functions were compared with that for bulk liquid water, and with a total radial distribution function determined experimentally for the two-layer hydrate of Na-montmorillonite by Powell et al. (1997). This comparison suggested that water molecules in the two-layer hydrate of montmorillonite have nearest-neighbor configurations which differ significantly from the tetrahedral ordering of nearest neighbors that characterizes bulk liquid water. Key Words--Adsorbed Water, Montmorillonite, Neutron Diffraction, Water. I N T R O D U C T I O N Neutron diffraction studies of interlayer water in hydrated 2 : 1 clay minerals based on isotopic-difference techniques have provided valuable insight as to the molecular structure of the electrical double layer on the surfaces of these minerals (Hawkins and Egelstaff, 1980; Skipper et al., 1991, 1994, 1995c). Powell et al. (1997) recently improved on the pioneering work o f Hawkins and Egels taff (1980) with a H/D isotopicdifference invest igat ion of water in the two-layer hydrate of Na-Wyoming montmoril lonite . Their principal experimental result was a first-order-difference total radial distribution function for adsorbed D20 molecules, obtained under condit ions in which no D rep lacement of clay mineral H atoms occurred. This experimental distribution function was compared to a hypothetical model distribution function based on available neutron diffraction data for bulk l iquid water (Soper and Phillips, 1986). Powel l et al. (1997) noted significant deviations be tween the interlayer total radial distribution function for adsorbed water and their mode l over intermolecular distances be tween 1.5-3 A, leading them to conclude that the nearest-neighbor coordination o f adsorbed water molecules in the twolayer hydrate differed significantly f rom the well known local tetrahedral coordinat ion that exists among water molecules in the bulk l iquid (Beveridge et al., 1983; Kusal ik and Svishchev, 1994). This same conclusion was obtained in a variety of recent Monte Carlo simulation studies of the molecular structure of interlayer water in the two-layer hydrate of Na -Wyoming montmori l loni te (Chang et al. , 1995; Boek et al. , 1995a, 1995b; Karaborni e t al. , 1996). Chang et al. (1995) and Boek et al. (1995b) both noted that water molecules in this hydrate were organized primari ly to accommodate the solvation requirements of Na + counterions bound in surface complexes on the basal planes of the clay mineral. Karaborni et al. (1996) commented that water protons were also attracted to surface oxygen ions near tetrahedral charge sites on the clay mineral to form hydrogen bonds, whereas water oxygens were attracted to the protons in structural OH groups at the bot tom of the ditr igonal cavities in the clay mineral surface. This compet i t ion for adsorbed water molecules by counterions, surface oxygen ions, and OH protons in the clay mineral structure would tend to disrupt the local tetrahedral coordination among water molecules that characterizes the bulk liquid. Chang et al. (1995) calculated radial distribution functions representing O-O and O-H spatial correlations in interlayer water on Na-montmori l loni te as part of the interpretation of their Monte Carlo simulation output. These two radial distribution functions, along with that for H-H spatial correlations, are in fact the principal contributors to the experimental total radial distribution function reported by Powel l et al. (1997). The present study was undertaken, therefore, to perform Monte Carlo simulations of the total radial distribution function for interlayer water in the two-layer hydrate of Na-montmori l loni te . The total radial distribution function for interlayer water in K-montmori l lonite also was simulated to evaluate more precisely the effects of cation solvation, which is much weaker for K + than it is for Na + (Ohtaki and Radnai, 1993). To examine the inherent quality of our simulations of interlayer water structure in more detail, bulk l iquid O-O, O-H, and H-H radial distribution functions based on the mode l for water molecule interactions used in our Monte Carlo simulations were compared to recent Copyright 9 1999, The Clay Minerals Society 192 Vol. 47, No. 2, 1999 Monte Carlo simulation of interlayer water in montmorillonite 193 experimental O-O, O-H, and H-H radial distribution functions for liquid water based on H/D isotopic-difference neutron diffraction data (Soper et al., 1997). TOTAL RADIAL DISTRIBUTION FUNCTION The measured intensity of neutron diffraction by an aqueous system, after instrumental calibration and corrections for unwanted scattering events, yields information about the spatial distribution of atoms in the system through the relationship (Enderby and Neilson, 1981): N (R sin kr Fr(k) = 4~r-~ Jo G(r)~-7--r r2 dr (1) where NFT(k) is the intensity of neutron diffraction at wavenumber k > 0, N is the total number of atoms in a system with volume V -= 4"rrR3/3 and G(r) = ~ ~ c,c~b~b~[g,~(r) 1] (2) is termed the total radial distribution function (Soper, personal communication, 1998), whose Fourier transform in Equation (1) is proportional to Fr(k). The experimental context for Equation (1), isotopic-difference neutron diffraction, is reviewed by Skipper et al. (1991, 1994) for applications to hydrated 2 :1 clay minerals, and by Enderby (1983), Neilson and Enderby (1989), Skipper and Nielsen (1989) and Ohtaki and Radnai (1993) for applications to aqueous solutions. The total radial distribution function G(r) provides a relative measure of the likelihood of coherent neutron scattering at wavenumber k by atoms located at a distance r from some atom placed at the origin of coordinates (Enderby and Neilson, 1981). Equation (2) shows that G(r) depends on the atomic fraction of the scattering atom c~, on its coherent scattering length b~ and on the radial distribution function g~(r), which is defined implicitly by the equation (Enderby and Neilson, 1981): dn,~ = 41TN~g~o(r)r 2 dr. (3) V " In Equation (3), dn~ is the average number of B-species atoms within a spherical shell of radius r and thickness dr enclosing an a-species atom which has been placed at r = 0, and N~ is the total number of [3 atoms in the system of volume V. Thus g~(r) is the relative probability that a [3 atom resides within dr at a radial distance r from an c~ atom centered at the origin of coordinates (Allen and Tildesley, 1987). It is normalized such that N~(p) = 4~r V g,~(r)r 2 dr (4) is the number of 13 atoms around the central c~ atom within a sphere of radius O ~ R (Beveridge et aL, 1983). As r "~ R, g~(r) -1 in an isotropic system (where $ means "approaches from below"), signifying a uniform distribution of [3 atoms around the central a atom, and G(r) in Equation (2) correspondingly goes to zero, indicating no likelihood of coherent neutron scattering by atoms far from that placed at the origin of coordinates. The physical role played by G(r) in Equation (2) is as a relative probability for finding any 13 atom that appears non-randomly at a distance r from an a atom placed at the origin, depending on the relative population of each type of atom (c~ or c~) and on its ability to scatter neutrons coherently (b~ or b~). Effective contributions to G(r) and, therefore, to neutron diffraction, come from atoms that are abundant, with large coherent scattering lengths, and which have strong spatial correlations with atoms of any of the species considered. In applications to neutron diffraction by the atoms in the interlayer region of a hydrated clay mineral, Equations (1) and (2) implicitly must be averages over spatial correlations that depend on the direction from a central atom, as well as on the intermolecular distance, given the constrained geometry of clay interlayers. Thus, G(r) and its component radial distribution functions are interpreted structurally as angularly-averaged quantities when applied to species in clay interlayers.